---
title: "Schäfer_BJPS_figures & tables"
author: "Armin Schäfer, WWU Münster"
date: |
  `r lubridate::today()`
output: html_document
---

This R Markdown (Rmd) file contains material to reproduce the figures and tables in Schäfer **"Cultural Backlash? How (not) to explain the rise of authoritarian populism."**. Please, run the files "prep_ess1_7.Rmd", "prep_ess9.Rmd" and "prep_bes.Rmd" first and save the created data in an appropriate folder.

```{r setup, include=FALSE}

knitr::opts_chunk$set(echo = TRUE, warning = FALSE, message = FALSE)

```

# Load packages and data----

*Note*: I rely on the [R package management tool `pacman`](). Please run `install.packages("pacman")` to install the package from CRAN. The command `pacman::p_load()` will then import all specified packages and install them first in case they are currently not installed on your machine.

```{r load packages and data, message = FALSE}
# Loading R-libraries.

if (!require("pacman")) install.packages("pacman")
pacman::p_load(tidyverse, dplyr, naniar, here, scales, 
               corpcor, modelsummary, ggdist, patchwork,
               ggpubr, readxl, flextable, officer, report)

# useful function
tabl <- function(...) table(..., useNA="ifany")

# load the datasets created with "

ess_1_7 <- readRDS(here("processed-data", "ess_1_7.rds"))
ess9 <- readRDS(here("processed-data", "ess9.rds"))
bes <- readRDS(here("processed-data", "bes.RDS")) 

```

## Median age per cohort
```{r median age}

ess_1_7 %>%
  group_by(cohort)%>%
  summarise(medage = median(age, na.rm=T))

```

# Figures in the paper----

## Figure_1

Example: Gays and lesbians should be free to life their lives at the wish
```{r fig_1, message = FALSE, warning = FALSE}

p1 <- ess_1_7 %>%
  filter(across(one_of(c("weight", "freehms", "cohort")),
                ~ !is.na(.x))) %>%
  mutate(freehmsz = wt.scale(freehms, weight)) %>%
  group_by(cohort) %>%
  summarise(mfree = weighted.mean(freehmsz, weight, na.rm=TRUE)) %>%
ggplot()+
  geom_hline(yintercept=0, color="grey70")+
  geom_point(aes(x = cohort, y = mfree),
             size=3)+
  labs(x = NULL, y = "Z-transformed scale\n")+
  scale_x_discrete(breaks=c("1", "2", "3", "4"),
                   labels=c("Interwar", "Boomers", 
                            "Gen X", "Millennials"))+
  ylim(-0.5,0.5)+
  theme_bw()+
  theme(axis.text = element_text(color="black", size=10))


p2 <- ess_1_7 %>%
  filter(!is.na(cohort), !is.na(weight)) %>%
  group_by(cohort) %>%
  summarise(mfree = weighted.mean(freehms, weight, na.rm=TRUE)) %>%
ggplot()+
  geom_point(aes(x = cohort, 
                 y = mfree),
             size=3)+
  geom_hline(yintercept=3, color="grey70")+
  ylim(1,5)+
  labs(x = NULL, y = "Original scale\n")+
  scale_x_discrete(breaks=c("1", "2", "3", "4"),
        labels=c("Interwar", "Boomers", 
                 "Gen X", "Millennials"))+
  theme_bw()+
  theme(axis.text = element_text(color="black", size=10))

p3 <- ggarrange(p1, p2)

p3

```


## Figure_2

Authoritarian and libertarian values---z-transformed factors
```{r fig_2}

p1_1 <- ess_1_7 %>%
  filter(across(one_of(c("weight", "RC1", "RC2", "cohort")),
                ~ !is.na(.x))) %>%
  mutate(authz = wt.scale(RC1, weight),
         libz = wt.scale(RC2, weight)) %>%
  group_by(cohort) %>%
  summarise(authoritarian = weighted.mean(authz, weight, na.rm=T),
            libertarian  = weighted.mean(libz, weight, na.rm =T)) %>% 
  select(cohort, authoritarian, libertarian) %>%
  pivot_longer(-cohort, names_to = "scale", values_to = "values") %>%
ggplot()+
  geom_hline(yintercept=0,
             color = "grey70")+
  geom_line(aes(x = cohort, y = values, 
                group=scale,
                linetype=scale), size=1)+
  geom_point(aes(x = cohort, y = values), size=2)+
  labs(x = NULL, y = "z scores",
       group = NULL, linetype = NULL)+
  scale_x_discrete(breaks=c("1", "2", "3", "4"),
        labels=c("Interwar", "Boomers", "Gen X", "Millennials"))+
  scale_linetype_manual(values=c("twodash", "dotted"))+
  ylim(-.5,.5)+
  theme_bw()+
  theme(axis.text = element_text(color="black",
                                 size = 10))
p1_1

```

## Figure_3

Authoritarian and libertarian values---original scales.
```{r fig_3}

p3_1 <- ess_1_7 %>%
  select(cohort, weight, behave, secure, safety, tradition, rules) %>%
  filter(across(one_of(c("weight", "cohort")),
                ~ !is.na(.x))) %>%
  pivot_longer(-(1:2), names_to="item", values_to="index") %>%
  group_by(cohort, item) %>%
  summarise(mindex = weighted.mean(index, weight, na.rm=T)) %>% 
ggplot()+
  geom_hline(yintercept = 3.5, color = "grey70")+
  geom_point(aes(x = cohort, y = mindex))+
  scale_x_discrete(breaks=c("1", "2", "3", "4"),
                   labels=c("Interwar", "Boomers", 
                            "Gen X", "Millennials"))+
  ylim(1,6)+
  labs(x = NULL,
       y = "Authoritarian values")+
  theme_bw()+
  theme(axis.text = element_text(color="black", size=8))+
  coord_flip()+
  facet_wrap(~item, ncol=5)


p3_2 <- ess_1_7 %>%
  filter(across(one_of(c("weight", "cohort")),
                ~ !is.na(.x))) %>%
  select(cohort, weight, surprise, adventure, creative, free, listen) %>%
  pivot_longer(-(1:2), names_to="item", values_to="index") %>%
  group_by(cohort, item) %>%
  summarise(mindex = weighted.mean(index, weight, na.rm=TRUE)) %>% 
ggplot(aes(x = cohort, y = mindex))+
  geom_hline(yintercept = 3.5, color="grey70")+
  geom_point(aes(x = cohort, y = mindex))+
  scale_x_discrete(breaks=c("1", "2", "3", "4"),
                   labels=c("Interwar", "Boomers", "Gen X", "Millennials"))+
  ylim(1,6)+
  labs(x= NULL,
       y = "Libertarian values")+
  theme_bw()+
  theme(axis.text = element_text(color="black", size=8))+
  coord_flip()+
  facet_wrap(~item, ncol=5)

p3_3 <- ggarrange(p3_1, p3_2, nrow=2)

p3_3

```

## Figure_4

Populist values (sum of trstplt, trstprl, trstprt re-scaled to 100).
```{r fig_4}

by_cohort <- ess_1_7 %>%
  filter(across(one_of(c("weight", "cohort")),
                ~ !is.na(.x))) %>%
  group_by(cohort) %>%
  summarise(mpop = weighted.mean(populist, weight, na.rm=TRUE))

plotpop1 <- ess_1_7 %>%
  filter(across(one_of(c("weight", "cohort")),
                ~ !is.na(.x))) %>%
  group_by(cohort, cntry) %>%
  summarise(mpop = weighted.mean(populist, weight, na.rm=TRUE)) %>%
ggplot()+
  geom_jitter(aes(x = cohort, y = mpop,
                  color=cntry),
              color="grey50", 
             size=2,
             width=0.1,
             alpha = 0.5)+
  geom_point(data = by_cohort, aes(x=cohort, y=mpop),
             size=3.5,
             color="black")+
  labs(x = NULL, y = "Populist attitudes\n (reversed political trust scale)")+
  geom_hline(yintercept=50, color="grey70")+
  ylim(0,100)+
  scale_x_discrete(breaks=c("1", "2", "3", "4"),
                   labels=c("Interwar", "Boomers",
                            "Gen X", "Millennials"))+
  theme_bw()+
  theme(axis.text = element_text(color="black", size=10))

plotpop1

```

## Figure_5

Authoritarian vs. libertarian values and populist attitudes in Great Britain. Data: British Election Study 2014-2023 Internet Panel, wave 9.
```{r fig_5}

half1 <- bes %>%
  filter(!is.na(cohort)) %>%
  filter(!is.na(al_scaleW7_W9)) %>%
  ggplot(aes(y = cohort, x = al_scaleW7_W9)) +
  stat_halfeye()+
  labs(y = NULL,
       x = "Libertarian-authoritarian scale")+
  theme_minimal()+
  theme(axis.text = element_text(color = "black"),
        legend.position="")


half2 <- bes %>%
  filter(!is.na(cohort)) %>%
  filter(!is.na(populism2)) %>%
  ggplot(aes(y = cohort, x = populism2)) +
  stat_halfeye()+
  labs(y = NULL,
       x = "Populism scale")+
  theme_minimal()+
  theme(axis.text = element_text(color = "black"),
        legend.position = "")


half1 / half2

```

# Tables in the paper----

## Table 1
```{r table-1, warning=FALSE, message = FALSE}

models <- list()
models[['Model 1']] <- lm(populist~ cohort + schwartzauth + factor(cntry), data = ess_1_7)
models[['Model 2']] <- lm(populist~ cohort + schwartzauth + factor(cntry) +
                            factor(essround), data = ess_1_7)
models[['Model 3']] <- lm(populist~ cohort + schwartzauth + factor(cntry) +
                            gender + factor(class5) + domicil + relig +
                            factor(essround), data = ess_1_7)

cm <- c('cohort2' = 'Baby Boomers',
        'cohort3' = 'Generation X',
        'cohort4' = 'Millennials',
        'schwartzauth' = 'Authoritarianism',
        'gender' = 'Gender',
        'factor(class5)2' = 'Lower middle class',
        'factor(class5)3' = 'Small business owners',
        'factor(class5)4' = 'Skilled workers',
        'factor(class5)5' = 'Low-skilled workers',
        'domicil' = 'Domicile',
        'relig' = 'Religiosity',
        '(Intercept)' = 'Constant')

rows <- tribble(~term, ~Model1, ~Model2, ~Model3,
                'Interwar generation (ref)', '',   '', '',
                'Upper middle class (ref)', '',   '', '',
                'Country fixed effects', 'yes', 'yes', 'yes',
                'Year fixed effects', 'no', 'yes', 'yes')

attr(rows, 'position') <- c(1,12,27,28)

modelsummary(models, output = "html", 
         coef_omit = 'cntry|essround',
         title ="Regression coefficients for low political trust (proxy for populism)",
         stars = TRUE,
         vcov = ~ cntry,
         coef_map = cm,
         add_rows = rows,
         gof_omit = 'DF|Deviance|AIC|BIC|Adj.R2|Log.Lik.',
         notes = list('Data: European Social Survey, rounds 1 to 7.',
                      'Note: The table shows the coefficients of OLS regression models that estimate populist attiudes (=low trust in political institutions).'))

```

## Table 2
```{r table-2, warning=FALSE, message = FALSE}

models <- list()
models[['Authoritarian']] <- lm(auth~ cohort + factor(cntry), weights=weight, data = ess9)
models[['Populist']] <- lm(pop~ cohort, weights=weight, data = ess9)
models[['Auth. Pop. (continuous)']] <- lm(popauth~ cohort + factor(cntry), weights=weight, data = ess9)
models[['Auth. Pop. (binary)']] <- glm(authpop~ cohort + factor(cntry), weights=weight, data = ess9, family=binomial())


rows2 <- tribble(~term, ~m1, ~m2, ~m3, ~m4,
                'Country fixed effects', 'yes', 'yes', 'yes', 'yes')

attr(rows2, 'position') <- c(9)

cm <- c('cohortBaby Boomers' = 'Baby Boomers',
        'cohortGen X' = 'Generation X',
        'cohortMillennials' = 'Millennials',
        '(Intercept)' = 'Constant')

msummary(models, output = 'html', 
         stars = TRUE,
         title = "Regression coefficients for voting behavior, ESS 9",
         vcov = ~ cntry,
         coef_map = cm,
         gof_omit = 'DF|Deviance|AIC|BIC|Adj.R2|Log.Lik.',
         add_rows = rows2,
         notes = list('Data: European Social Survey, round 9', 
         'Note: The table shows the coefficients of four regression models that estimate the vote for authoritarian, populist and authoritarian populist parties.')) 

```

# Appendix

## Table A-1
```{r appendix-A-1, warning=FALSE, message = FALSE}

variables <- read_excel(here("processed-data","ess_variable_coding.xlsx"))

flextable(variables) %>%
  width(j = 1, width = 1.5) %>%
  width(j = 2, width = 3) %>%
  width(j = 3, width = 2) %>%
  line_spacing(space = 1, part = "all") %>%
  hline_top(border = fp_border(width=1.1, color ="black"), part="header") %>%
  hline(border = fp_border(width = .5, color = "grey50"), part = "body" ) %>%
  hline_bottom(border = fp_border(width=1.1, color ="black"), part="body") %>%
  hline(border = fp_border(width = 1.1, color = "black"), part = "header") %>%
  set_caption("Coding of variables, European Social Survey")

```

## Table A-2
```{r appendix-A-2, echo=FALSE, warning=FALSE, ft.align="left"}

ess1_7 <- ess_1_7 %>% 
  select(cohort, domicil, relig, populist, schwartzauth,
         behave, rules, safety, secure, 
         tradition, adventure, creative,
         free, listen, surprise,
         imueclt, euftf, imwbcnt) %>%
  rename(., 
         Cohort = cohort,
         Domicile = domicil,
         Religiosity = relig,
         "Authoritarian values" = schwartzauth,
         "Populist (low trust)" = populist,
         Behave = behave,
         Rules = rules,
         Safety = safety,
         Secure = secure,
         Tradition = tradition,
         Adventure = adventure,
         Creative = creative,
         Free = free,
         Listen = listen,
         Surprise = surprise,
         "EU integration" = euftf,
         "Mig. culture" = imueclt,
         "Mig. economy" = imwbcnt) %>%
  mutate_all(as.numeric)

app1 <- datasummary(All(ess1_7) ~ Mean + SD + Min + Max + N,
                    data = ess1_7,
                    notes = c('Data: European Social Survey 1-7.'),
                    output = 'flextable') 


app1 %>%
  width(width = 1) %>%
  width(j = 1, width = 2.5) %>%
  line_spacing(space = 1, part = "all") %>%
  hline_top(border = fp_border(width=1.1, color ="black"), part="header") %>%
  hline(border = fp_border(width = .5, color = "grey50"), part = "body" ) %>%
  hline_bottom(border = fp_border(width=1.1, color ="black"), part="body") %>%
  hline(border = fp_border(width = 1.1, color = "black"), part = "header") %>%
  set_caption("Descriptive statistics, EES 1-7")

```

## Table A-3

```{r appendix-A-3, echo=FALSE, warning=FALSE, ft.align="left"}

ess_select <- ess9 %>% 
  select(gender, cohort, married, separate, children,
         edu5, ethnic, relig, diseco,
         schwartzauth, trust, efficacy, immig) %>%
  rename(., 
         Gender = gender,
         Cohort = cohort,
         Married = married,
         Separated = separate,
         Children = children,
         Education = edu5,
         "Ethnic minority" = ethnic,
         "Satisfaction economy" = diseco,
         Religiosity = relig,
         "Authoritarian values" = schwartzauth,
         "Populist values" = trust,
         "External efficacy" = efficacy,
         "Pro immigration" = immig)


app2 <- datasummary(All(ess_select) ~ Mean * Arguments(fmt = "%.2f") + 
            SD  * Arguments(fmt = "%.2f") + 
            Min + Max + (N = 1),
            fmt = "%.0f",
            data = ess_select,
            notes = c('Data: European Social Survey 9.'),
            output = 'flextable')

app2 %>%
  width(width = 1) %>%
  width(j = 1, width = 2.5) %>%
  line_spacing(space = 1, part = "all") %>%
  hline_top(border = fp_border(width=1.1, color ="black"), part="header") %>%
  hline(border = fp_border(width = .5, color = "grey50"), part = "body" ) %>%
  hline_bottom(border = fp_border(width=1.1, color ="black"), part="body") %>%
  hline(border = fp_border(width = 1.1, color = "black"), part = "header") %>%
  set_caption("Descriptive statistics, EES 9")


```

## Table A-4

```{r appendix-A-4, echo=FALSE, warning=FALSE, ft.align="left"}

ess9_factor <- ess9 %>%
  select(cohort, class5, unemp) %>%
    rename(., 
           "Ever unemployed" = unemp,
           "Class scheme" = class5,
           "Birth cohort" = cohort)


app3 <- datasummary_skim(
  ess9_factor,
  type = "categorical",
  fmt = "%.1f",
  notes = c('Data: European Social Survey 9.'),
            output = 'flextable')

app3 %>%
  width(j = 2, width = 4) %>%
  #width(width = 1) %>%
  line_spacing(space = 1, part = "all") %>%
  hline_top(border = fp_border(width=1.1, color ="black"), part="header") %>%
  hline(border = fp_border(width = .5, color = "grey50"), part = "body" ) %>%
  hline_bottom(border = fp_border(width=1.1, color ="black"), part="body") %>%
  hline(border = fp_border(width = 1.1, color = "black"), part = "header") %>%
  set_caption("Descriptive statistics, EES 9, factors")

```

## Table A-5

```{r appendix-A-5, echo=FALSE, warning=FALSE, message = FALSE, results="asis"}

models <- list()
models[['Model 1']] <- lm(popauth~ cohort + factor(cntry), weights=weight, data = ess9)

models[['Model 2']] <- lm(popauth~ cohort + gender + married + separate + children + 
                            class5 + edu5 + ethnic + relig + unemp + lrscale +
                            diseco + factor(cntry), weights=weight, data = ess9)

models[['Model 3']] <- lm(popauth~ cohort + gender + married + separate + children + 
                            class5 + edu5 + ethnic + relig + unemp + lrscale +
                            diseco + schwartz01 + trust01 + immig01 + 
                            factor(cntry), weights=weight, data = ess9)

models[['Model 4']] <- lm(popauth~ cohort + gender + married + separate + 
                            children + 
                            class5 + edu5 + ethnic + relig + unemp + lrscale +
                            diseco + schwartz01 + efficacy01 + immig01 + 
                            factor(cntry), weights=weight, data = ess9)

row1 <- c('Country fixed effects', 'yes', 'yes', 'yes', 'yes')

cm <- c('cohortBaby Boomers' = 'Baby Boomers',
        'cohortGen X' = 'Generation X',
        'cohortMillennials' = 'Millennials',
        'children' = 'Children',
        'married' = 'Married',
        'separate' = 'Separate',
        'gender' = 'Gender',
        'class5Unskilled workers' = 'Unskilled workers',
        'class5Skilled workers' = 'Skilled workers',
        'class5Lower-grade service class' = 'Lower-grade service class',
        'class5Small business owners' = 'Small business owners',
        'edu5' = 'Education',
        'ethnic' = 'Ethnic minority',
        'domicil' = 'Domicile',
        'relig' = 'Religiosity',
        'uemp12m' = 'Ever been unemployed for 12 months',
        'hincfel' = 'Subjective financial insecurity',
        'diseco' = 'Dissatisfied with state of economy',
        'lrscale' = 'Placement on left-right scale',
        'schwartz01' = 'Authoritarianism',
        'efficacy01' = 'External efficacy',
        'trust01' = 'Political trust',
        'immig01' = 'Pro immigration',
        '(Intercept)' = 'Constant')

attr(rows2, 'position') <- c(43)

table3 <- msummary(models, output = 'flextable', 
         title = "Regression coefficients: vote for authoritarian-populist parties, ESS9",
         stars = TRUE,
         vcov = ~ cntry,
         coef_map = cm,
         add_rows = rows2,
         gof_omit = 'DF|Deviance|AIC|BIC|Adj.R2|Log.Lik.',
         notes = list('Note: The table show the coefficients of OLS regression models that estimate the vote for authoritarian-populist parties.'))

table3 %>%
  width(width = 1.2)%>%
  width(j = 1, width = 2)%>%
  line_spacing(space = 0.7, part = "all") %>%
  hline_top(border = fp_border(width=1.1, color ="black"), part="header") %>%
  hline(border = fp_border(width = .5, color = "grey50"), part = "body" ) %>%
  hline_bottom(border = fp_border(width=1.1, color ="black"), part="body") %>%
  hline(border = fp_border(width = 1.1, color = "black"), part = "header")


```

## Figure A-1

Attitudes in comparison: different scales
```{r fig_A-1}

p2_1 <- ess_1_7 %>%
  mutate(euftf = as.numeric(euftf)) %>%
  filter(across(one_of(c("cohort", "imueclt", "imwbcnt", "euftf", "weight")),
                ~ !is.na(.x))) %>%
  select(cohort, imueclt, imwbcnt, euftf, weight)%>%
  mutate(imuecltz = wt.scale(imueclt, weight),
         euftfz = wt.scale(euftf, weight),
         imwbcntz = wt.scale(imwbcnt, weight)) %>%
  group_by(cohort) %>%
  summarise("Mig. culture" = weighted.mean(imuecltz, weight, na.rm=TRUE),
            "EU integration" = weighted.mean(euftfz, weight, na.rm=TRUE),
            "Mig. economy" = weighted.mean(imwbcntz, weight, na.rm=TRUE)) %>%
  pivot_longer(-cohort, names_to = "items", values_to = "index") %>%
ggplot()+
  geom_hline(yintercept = 0, color="grey70")+
  geom_point(aes(x = cohort, y = index), size=2)+
  labs(y = "z-transformed scores\n",
       x = NULL,
       color = NULL)+
  scale_x_discrete(breaks=c("1", "2", "3", "4"),
                   labels=c("Interwar", "Boomers", 
                            "Gen X", "Millennials"))+
  ylim(-0.3,0.3)+
  theme_bw()+
  theme(axis.text = element_text(color = "black", size=8),
        legend.position = "top")+
  facet_wrap(~items)


p2_2 <- ess_1_7 %>%
  mutate(euftf = as.numeric(euftf)) %>%
  filter(across(one_of(c("cohort", "imueclt", "imwbcnt", "euftf", "weight")),
                ~ !is.na(.x))) %>%
  select(cohort, imueclt, imwbcnt, euftf)%>%
  group_by(cohort) %>%
  summarise("Mig. culture" = weighted.mean(imueclt, na.rm=TRUE),
            "EU integration" = weighted.mean(euftf, na.rm=TRUE),
            "Mig. economy" = weighted.mean(imwbcnt, na.rm=TRUE)) %>%
  pivot_longer(-cohort, names_to = "items", values_to = "index") %>%
ggplot()+
  geom_hline(yintercept=5, color="grey70")+
  geom_point(aes(x = cohort, y = index), size=2)+
  labs(y = "Original scale\n",
       x = NULL,
       color = NULL)+
  scale_x_discrete(breaks=c("1", "2", "3", "4"),
                   labels=c("Interwar", "Boomers", 
                            "Gen X", "Millennials"))+
  ylim(0,10)+
  theme_bw()+
  theme(axis.text = element_text(color = "black", size=8),
        legend.position = "top")+
  facet_wrap(~items)

p2_3 <- ggarrange(p2_1, p2_2, nrow=2)

p2_3

```

## Packages used

```{r packages}
cite_packages()

```

